A method of estimating blood pressure of a subject

ABSTRACT

The present invention relates to a method and system for estimating blood pressure of a subject. In particular, but not exclusively, the method involves receiving a photoplethysmogram (PPG) signal from a light-based Pulse-Plethysmography sensor applied to the skin of a subject and reconstructing a pulse blood pressure waveform between systolic and diastolic blood pressure of the subject. Additionally, but not exclusively, the method involves processing the pulse blood pressure waveform and reconstructing an absolute blood pressure waveform of the subject.

TECHNICAL FIELD

The present invention relates to a method and system for estimatingblood pressure of a subject. In particular, but not exclusively, themethod involves receiving a photoplethysmogram (PPG) signal from alight-based Pulse-Plethysmography sensor applied to the skin of asubject and reconstructing a pulse blood pressure waveform betweensystolic and diastolic blood pressure of the subject. Additionally, butnot exclusively, the method involves processing the pulse blood pressurewaveform and reconstructing an absolute blood pressure waveform of thesubject.

BACKGROUND OF INVENTION

Blood pressure (BP) measurement in health care settings and day-to-daylife helps in monitoring and diagnosing a multitude of physiological andpathological responses. Current non-invasive BP measurement systems havetypically been focussed on obtaining the systolic (SBP) and diastolicblood pressure (DBP) values of a subject. These systems typically usedigital oscillometry techniques or various direct equation-basedapproaches attained from bio-signals using, for example,electrocardiogram (ECG), photoplethysmogram (PPG), pressure andsound-based sensors.

Oscillometry techniques are commonly used with pressure sensors. Forexample, a pump is used to inflate a rubber cuff on a subject and theresults of the pressure sensor are then used to estimate the SBP and DBPlevels of the subject. These techniques, however, are inaccurate andresults may vary up to ±30 mmHg from the real pressure value. Thisinaccuracy can largely be attributed to the variance in the underlyingalgorithms used, which may not consider dynamic physiological changes.Further, devices employing these techniques are bulky and infeasible forportable usability. Also, 24-hour monitoring is hindered due to thediscomfort of their continued use on the subject as a result of thefrequent compression of the cuff required around the desired pressuremeasurement location.

Equation-based approaches have been under research for a significantamount of time. Non-invasive BP measurement using these approaches hastypically been implemented in research studies rather than in commercialuse. These approaches use various sensing systems to acquire bio-signalsfrom a subject to develop a polynomial or logarithmic relationship withthe SBP and DBP. The accuracy of these approaches, however, hastypically been compromised due to the inefficiency of the equations inconsidering the many factors that affect blood pressure. Further,frequent calibration of devices employing these approaches is typicallyrequired to maintain an accurate measurement over time.

In particular, standard mathematical regression models typically fail tocreate a uniform relationship across all subjects. Also, in research,the sample size selected is typically inadequate, leading to biasedmodels. The focus might be limited to a control group which will notgeneralise for an entire population including healthy subjects.Mathematical models may also focus on the present sensor estimates anddo not consider the short term or the long-term physiological responsesof the subjects, such as a subject suffering from an elevated pressuredue to the presence of a medical professional. This may lead to falsealarms or false diagnosis.

The amount of data required to create generalised mathematical models isobviously very large and it requires a very large amount ofcomputational power compared to more specialised models. Such largedatabases may not be present or may be hard to access due to ethicalreasons. Furthermore, a researcher is typically required to hand-selectfeatures from sensor signals based on the researcher’s knowledge andfindings for these models. But, these features may or may not have anycorrelation. Moreover, it is difficult to identify features from thesensor signals correlating to mean arterial pressure (MAP) of a subjectdue to, for instance, physiological differences in test subjects.Further, missing features may further lead to inaccurate estimations.

Accordingly, the current gold-standard in obtaining BP measurements isto take invasive measurements of the BP of a subject using a catheter.Invasive measurements are, however, undesirable to the subject forobvious reasons and, in some cases, may still fail to meet the requiredlevel of accuracy of within 5 ±8 mmHg set by Food and DrugAdministration (FDA).

SUMMARY OF INVENTION

According to one aspect of the present invention, there is provided amethod of estimating blood pressure of a subject, the method including:receiving a photoplethysmogram (PPG) signal from a light-basedPulse-Plethysmography sensor applied to the skin of the subject;processing the PPG signal using a Wavelet transform algorithm to derivePPG wavelet coefficients in a plurality of frequency bands; estimatingblood pressure coefficients by processing the PPG wavelet coefficientsusing a machine learning algorithm that has been trained on trainingdata including PPG wavelet coefficients derived from PPGs fromlight-based Pulse-Plethysmography sensors applied to the skin of testsubjects correlated with Invasive Arterial Blood Pressure (ABP) waveletcoefficients derived from simultaneously obtained ABP measurements ofsystolic and diastolic blood pressure of the test subjects; andreconstructing a pulse blood pressure waveform between systolic anddiastolic blood pressure of the subject from the blood pressurecoefficients.

According to another aspect of the present invention, there is provideda method of monitoring blood pressure of a subject, the methodincluding: applying a light-based Plethysmography sensor to the skin ofthe subject for a period of time; and estimating the blood pressure ofthe subject according to the above method at designated intervals overthe period time.

According to another aspect of the present invention, there is provideda system for estimating blood pressure of a subject, the systemincluding: a light-based Pulse-Plethysmography sensor configured to beapplied to the skin of the subject to generate a photoplethysmogram(PPG) signal; a processor in data communication with the sensor; amemory; and software resident in the memory and accessible to theprocessor, the software including a series of instructions executable bythe processor to configure the processor to: process the PPG signalusing a Wavelet transform algorithm to derive PPG wavelet coefficientsin a plurality of frequency bands; estimate blood pressure coefficientsby processing the PPG wavelet coefficients using a machine learningalgorithm that has been trained on training data including PPG waveletcoefficients derived from PPGs from light-based Pulse-Plethysmographysensors applied to the skin of test subjects correlated with InvasiveArterial Blood Pressure (ABP) wavelet coefficients derived fromsimultaneously obtained ABP measurements of systolic and diastolic bloodpressure of the test subjects; reconstruct a pulse blood pressurewaveform between systolic and diastolic blood pressure of the subjectfrom the blood pressure coefficients; and output the pulse bloodpressure waveform to a display.

The method uses a machine learning algorithm to estimate blood pressurecoefficients by processing derived PPG wavelet coefficients. Existingmachine learning techniques, such as standard regression models, supportvector machines, classification and regression trees (CART) and neuralnetworks, could be used to predict systolic and diastolic bloodpressures from a PPG signal. Doing so, however, would require animpractical and uncountable number of features for input to the machinelearning algorithms to achieve a high level of accuracy in predictingthe SBP and DBP. The justification for the use of these features maythus become random or heuristic. Further practical implementation ofsearching such features becomes computationally intensive and exhaustivetoo - missing features which may result in the incorrect estimation ofthe pressure values.

The power of machine learning can be used in extracting andcomprehending features automatically rather than manual extraction withhuman interference. From a physiological perspective, this is essentialdue to the quasi-static nature of the human body. Preferably, the methoduses a semi-supervised approach in blood pressure estimation, enablingan accurate estimation focussing the machine learning algorithm onaspects where human intelligence may overlook. The method isincorporated into a non-invasive, non-pressurised system for estimatingBP of a subject. This non-invasive, non-pressurised system enablesportable and continuous monitoring of the subject without hinderingdaily activity of the subject. Further, such a system can replace theabove described invasive systems in circumstances where, for example,monitoring of BP might be essential but not critical in certain areas ofhospitals.

In an embodiment, the light-based Pulse-Plethysmography sensor includesa light source having an emission wavelength and a photodiode having adetection wavelength, wherein the emission wavelength and the detectionwavelength are around an isosbestic wavelength where oxygenated anddeoxygenated blood absorbs same amount of light. For example, the lightsource is an Infrared (IR) light source. The reconstruction of the pulseblood pressure waveform is thus able to be achieved non-invasively.While this embodiment implements an isosbestic sensor, it will beappreciated by those persons skilled in the art that the method andsystem could be configured to use a sensor incorporating wavelengthsanywhere from the visible spectrum to the infra-red spectrum. Theisosbestic sensor at the infra-red wavelength is preferable as itcompensates for the equal absorption by oxygenated and deoxygenatedblood and is less susceptible to motion artefacts. An estimation of BPusing the isosbestic sensor may therefore be more confident about theunderlying physiological changes being homeostatic rather thaninterference from external motion.

In an embodiment, the light-based Pulse-Plethysmography sensor includestwo or more pairs of the light source and the photodiode, and the methodfurther includes determining a signal to noise ratio (SNR) of the PPGsignal from each of the pairs of the light source and the photodiode,and selecting one of the pairs of the light source and the photodiodefor providing the PPG signal with the best SNR.

In an embodiment, the method further includes receiving an inertialmeasurement signal from an inertial measurement unit (IMU) adjacent thelight-based Pulse-Plethysmography sensor, and compensating for motionartefacts when determining the SNR of the PPG signal from each of thepairs of the light source and the photodiode includes using the inertialmeasurement signal.

In an embodiment, the method further includes processing selected onesof the PPG wavelet coefficients using the machine learning algorithmbased on an energy level present at each of the frequency bandsexceeding a threshold level. The method may further include providingthe selected ones of the PPG wavelet coefficients in a 2-Dimensionalmatrix of features to the machine learning algorithm. Further, eachcolumn of the 2-Dimensional matrix is autoregressive incorporating nearpast and near future samples of the selected ones of the PPG waveletcoefficients. Alternatively, the method may provide the PPG waveletcoefficients as a 1-Dimensional vector incorporating the reconstructedwavelet of specific decomposed wavelet level of blood pressure. Further,multiple 1-Dimensional vectors of reconstructed wavelets of all desiredlevels of blood pressure in 0.5 - 8 Hz bandwidth may be provided.

In an embodiment, the machine learning algorithm is a Recurrent NeuralNetwork (RNN) algorithm, such as a Long short-term memory (LSTM) model.Preferably, the method further includes processing the PPG waveletcoefficients using multiple LSTM models.

That is, embodiments of the invention provide for the estimation andreconstruction of a pulse blood pressure waveform to be achievedalgorithmically using a Wavelet signal processing technique along with amemory-based machine learning algorithm of Long Short-Term Memory(LSTM). One benefit of using LSTM networks relates to the ability toaccount for dynamically changes in the body of the subject byconsidering short and long past states of the body under homeostasis.The abovementioned existing machine learning techniques or directequation-based approaches mostly depend on present values from thesensors and fail to incorporate dynamic changes. Further, the use of theWavelet transform algorithm involves obtaining components of the PPGsignal at different frequency bandwidths to develop a moretime-frequency based continuous model.

Further, embodiment of the invention using LSTM along with the Wavelettransform algorithm enable a time-frequency based analysis of the PPGsignal. The Wavelets algorithm decomposes the PPG signal into differentfrequency bands. The PPG wavelets coefficients in each of the frequencybands are organised into 2-Dimensional matrix representation which is afeature to the neural network model. LSTM neural network models aredeveloped to establish a correlation between the PPG waveletcoefficients in specific frequency bands with that of Invasive ArterialBlood Pressure (ABP) wavelet coefficients. These correlations are usedin time series reconstruction of ABP wavelet coefficients. In theembodiment, the method further includes reconstructing output of themultiple LSTM models to form the pulse blood pressure waveform using aninverse Wavelet transform algorithm.

In an embodiment, the Wavelet transform algorithm is a MaximallyOverlapped Discrete Wavelet Transform (MODWT) algorithm, and the methodfurther includes reconstructing the pulse blood pressure waveform usingan inverse MODWT algorithm. That is, embodiments of the method providesfor the calculation of a pulse pressure, ΔP, following sequential methodsteps involving signal processing techniques with MODWT and memory-basedrecurrent neural network (RNN) estimation.

In an embodiment, the method further includes pre-processing the PPGwavelet coefficients before processing using the machine learningalgorithm by lowpass filtering the PPG wavelet coefficients with acut-off frequency for the frequency bands. For example, the cut-offfrequency is 10 Hz.

In an embodiment, the method further includes: processing the pulseblood pressure waveform using a further Wavelet transform algorithm toderive pulse blood pressure wavelet coefficients in a plurality offrequency bands; extracting features from the pulse blood pressurewavelet coefficients in the plurality of frequency bands; estimatingmean arterial blood pressure (MAP), systolic blood pressure (SBP) ordiastolic blood pressure (DBP) coefficients by processing the featuresusing a further machine learning algorithm that has been trained on thetraining data; and reconstructing the MAP, SBP and or DBP waveforms fromthe MAP, SBP and or DBP coefficients, respectively.

Preferably, the method further includes combining the pulse bloodpressure waveform and one or more of the MAP, SBP or DBP waveforms toreconstruct an absolute blood pressure waveform of the subject.

In an embodiment, the further machine learning algorithm is aConvolutional Long short-term memory (ConvLSTM) network and the featuresare a 2-Dimensional feature matrix including temporal and spatialfeatures of the pulse blood pressure wavelet coefficients.

In an embodiment, the method further includes estimating intermediateMAP, SBP or DBP coefficients using the ConvLSTM network.

In an embodiment, the method further includes: selecting features of thePPG signal and the pulse blood pressure waveform; and estimating furtherintermediate MAP, SBP or DBP coefficients by processing the features ofthe PPG signal and the pulse blood pressure waveform using deep neuralnetworks, respectively, that have been trained on the training data.

In this embodiment, the method further includes: concatenating theintermediate and further intermediate MAP, SBP or DBP to form aconcatenated output vector; and passing the output vector throughmultiple layers of hidden neurons in an output network to generate theestimated MAP, SBP or DBP coefficients.

In an embodiment of the system, the processor is further configured to:process the pulse blood pressure waveform using a further Wavelettransform algorithm to derive pulse blood pressure wavelet coefficientsin a plurality of frequency bands; extract features from the pulse bloodpressure wavelet coefficients in the plurality of frequency bands;estimate mean arterial blood pressure (MAP), systolic blood pressure(SBP) or diastolic blood pressure (DBP) coefficients by processing thefeatures using a further machine learning algorithm that has beentrained on the training data; and reconstruct the MAP, SBP and or DBPwaveforms from the MAP, SBP and or DBP coefficients, respectively.Preferably, the processor is further configured to: combine the pulseblood pressure waveform and one or more of the MAP, SBP or DBP waveformsto reconstruct an absolute blood pressure waveform of the subject; andoutput the absolute blood pressure waveform to the display.

BRIEF DESCRIPTION OF DRAWINGS

Embodiments of the present invention will now be described withreference to the accompanying drawings, wherein:

FIG. 1 is a block diagram showing a system for estimating blood pressureof a subject in accordance with an embodiment of the present invention;

FIG. 2 is a flow chart showing a system for estimating blood pressure ofa subject in operation in accordance with an embodiment of the presentinvention;

FIG. 3 is a flow chart showing a method of training a machine learningalgorithm in accordance with an embodiment of the present invention;

FIG. 4 is a graph showing an energy distribution of MODWT coefficientsfor ABP;

FIG. 5 is an example of a Long short-term memory (LSTM) neuronstructure;

FIG. 6 is an example of a Deep LSTM neural network;

FIG. 7 is a flow chart showing a method of estimating blood pressure ofa subject in operation in accordance with an embodiment of the presentinvention;

FIG. 8 is a flow chart showing a method of estimating blood pressure ofa subject in operation in accordance with an embodiment of the presentinvention;

FIG. 9 is an example of a blood pressure waveform of a subject;

FIG. 10 is a block diagram showing a method of estimating blood pressureof a subject in accordance with an embodiment of the present invention;

FIG. 11 is a block diagram showing hardware for estimating bloodpressure of a subject in accordance with an embodiment of the presentinvention;

FIG. 12 is a block diagram showing a method of selecting a light sourceand a photodiode for estimating blood pressure of a subject inaccordance with an embodiment of the present invention;

FIG. 13 is a flow chart showing a method of estimating blood pressure ofa subject in operation in accordance with an embodiment of the presentinvention;

FIG. 14 is a flow chart showing a method of training a machine learningalgorithm in accordance with an embodiment of the present invention;

FIG. 15 is an example of a ConvLSTM neuron;

FIG. 16 is an example of a multi-layered ConvLSTM structure;

FIG. 17 is a flow chart showing a method of estimating blood pressure ofa subject in operation in accordance with an embodiment of the presentinvention;

FIG. 18 a is a diagram showing ConvLSTM sequencing for sequence 1;

FIG. 18 b is a diagram showing ConvLSTM sequencing for sequence 2;

FIG. 18 c is a diagram showing ConvLSTM sequencing for sequence N; and

FIG. 19 is a diagram of a Mean Arterial Pressure estimation of a patientover 10 minutes.

DETAILED DESCRIPTION

FIG. 1 shows an embodiment of a system 10 for estimating blood pressureof a subject non-invasively. The system 10 includes a light-basedPulse-Plethysmography sensor 12 that is configured to be applied to theskin of a subject to generate a photoplethysmogram (PPG) signal. Thesensor 12 is applied to the skin by being attached, using a suitableattachment means, to a desired location on the body, such as a finger orwrist, of the subject. The sensor 12 includes a light source and aphotodiode to detect light reflected from the subject. For example, thesensor 12 includes an Infrared (IR) LED light source with an emissionwavelength driven by an LED driver around an isosbestic wavelength whereoxygenated and deoxygenated blood absorbs the same amount of light. Thephotodiode also has a detection wavelength around the isosbesticwavelength.

The system 10 further includes a computer 14 in data communication withthe sensor 12. The data connection may be wired or wireless. Forexample, the PPG signal is received via Bluetooth or WIFI. The computer14 includes a processor 16 in data communication with a memory 18. Thememory 18 includes software that includes a series of instructionsexecutable by the processor 16 to configure the processor 16 to performa number of steps when the system 10 is in the implementation phase.

The processor 16 processes the PPG signal received wirelessly from thesensor 12 using a Wavelet transform algorithm, such as a MaximallyOverlapped Discrete Wavelet Transform (MODWT) algorithm, to derive PPGwavelet coefficients in a plurality of frequency bands. As mentioned,the PPG signal is acquired from the body of the subject at the range ofthe isosbestic wavelength. It is then pre-processed, re-sampled anddecomposed with the MODWT algorithm which is discussed below in moredetail.

The processor 16 then estimates blood pressure coefficients byprocessing the PPG wavelet coefficients using a machine learningalgorithm that has been trained on training data including PPG waveletcoefficients derived from PPGs from light-based Pulse-Plethysmographysensors applied to the skin of test subjects correlated with InvasiveArterial Blood Pressure (ABP) wavelet coefficients derived fromsimultaneously obtained ABP measurements of systolic and diastolic bloodpressure of the test subjects. The training of the machine learningalgorithm was conducted prior to the implementation phase and, in oneembodiment, is a neural network algorithm which gives a time-frequencyestimation of blood pressure coefficients at multiple frequencybandwidths.

The processor 16 then reconstructs a pulse blood pressure waveformbetween systolic and diastolic blood pressure of the subject from theblood pressure coefficients using an inverse Wavelet transformalgorithm, such as an inverse MODWT algorithm. Finally, the processor 16outputs the pulse blood pressure waveform to a display 20. The displayof the pulse blood pressure waveform on the display 20 enables the bloodpressure of the subject to be monitored by, for example, careprofessionals when the sensor 12 is applied to the skin of the subjectfor a period of time.

FIG. 2 shows another embodiment of the system 10 in operation in theimplementation phase. In this embodiment, the sensor 12 is shown toinclude a number of components that are used to generate the PPG signaland transmit it to the computer 14. That is, the senor 12 includes askin sensor interface for attaching the sensor 12 to the desiredlocation on the body. The sensor 12 includes an isosbestic LED, drivenby an LED driver, and an isosbestic photodiode with a photodiode signalreceiver which receives the PPG signal sampled at a specific samplingfrequency F_(s). The PPG signal is then transmitted using a wirelesstransmitter of the sensor 12 to the computer 14 to undertake theabovementioned processing including machine learning processing andreconstruction.

The computer 14 firstly undertakes pre-processing of the received PPGsignal, which involves low pass filtering with a cut-off frequency at 20Hz. In this frequency range, the expectation of higher physiologicalresponse is present. Higher frequency components in the signal are dueto the presence of external interference. A check may be furtherconducted if the signal is valid or not through a motion removalartefact algorithm. The result of the pre-processing step is anormalised PPG signal between 0 - 1.

The computer 14 then carries out a Maximally Overlapped Discrete WaveletTransform (MODWT) algorithm on the normalised PPG signal. This generatesdifferent levels of coefficients at different frequency bands. nspecific frequency bandwidths are selected which have maximuminformation about the physiological changes. The output vector from eachfrequency bandwidth is structured into a 2-Dimensional matrix. Thismatrix is designed to be autoregressive with an expectation of thecurrent value of prediction to be depended on ƒ past samples and ƒ - 1future samples, giving it a non-causal representation.

The restructured 2-D matrix at each frequency bandwidth is presented toa trained machine learning algorithm. In the embodiment, the machinelearning algorithm includes multiple Long short-term memory (LSTM)models. n LSTM models generate n vectors of wavelet coefficientsrepresenting an estimation of blood pressure wavelet coefficients. Thesecoefficients are then passed through an Inverse Maximally OverlapDiscrete Wavelet Transform (IMODWT) algorithm to generate areconstructed differential (pulse) pressure waveform, ΔP.

Training Phase

An embodiment of the training phase of the machine learning algorithm ofthe system 10 is shown in the flow chart 22 of FIG. 3 . As mentioned,the training data for the machine learning algorithm includes PPGwavelet coefficients derived from PPGs from test subjects correlatedwith Invasive Arterial Blood Pressure (ABP) wavelet coefficients derivedfrom simultaneously obtained ABP measurements. More specifically, theABP measurements and PPG signals recorded from several test subjectssimultaneously were obtained through an ethics approval process from theMIMIC-II database for the system 10. The de-identified waveforms fromthe test subjects were pre-processed and filtered. The initial samplingof the signals was conducted at 500 Hz. It was determined that such highsampling is not required as the desired interest of frequency bands werebelow 20 Hz. A resampling at 125 Hz was thus done for the signals. Itwas determined that the choice of the sampling frequency at 125 Hz hereis heuristic, and a lower sampling frequency can be chosen. Two signalswere structured, starting from the same timestamp.

The PPG signals were normalised as the recording from different patientswere undertaken using different devices having different voltage output.Cross-Correlation between the two signals was undertaken. Thecross-correlation ensured to estimate the phase lag between the twosignals.

r(t) = ABP(t) ⊗ PPG(t) = ∫_(−∞)^(∞)ABP(τ)PPG(τ + t)dτ

The cross-correlation estimate r(t) in equation (1) is used for thecalculation of the maximum lag,

ε = t[r(t)|_(max (r(t))))]

This phase lag, ε is generated because of the different locations of therecording undertaken for the pressure and the PPG in the patients. ThePPG signal is shifted forward (or backward) based on the lag. This wasessential as the lag is a variable depending on the distance between theplacement of the two sensors on the patient. Since this is an unknownvariable, it creates uncertainty in the estimates in blood pressure inmodel training.

MODWT Transform

Maximally overlap discrete wavelet transform (MODWT) is used for thedecomposition of the signal. Compared to the standard wavelet transform(DWT), MODWT undertakes a multi-resolution analysis of the signal. Thesize of the coefficients for all levels are defined for all samplesizes. No decimation occurs for subsequent frequency bands coefficients,unlike standard DWT.

The retainment of the same number of coefficients enables in developingan input-output time-series relationship in different frequencybandwidths. Whereas, DWT coefficients do not provide any time-seriescorrelation due to their time-varying property.

In MODWT, the signal is passed through a set of high pass and low passquadrature mirror filters. There have been several different types ofquadrature mirror filters proposed in theory such as Daubechies,Symlets, Coiflets, Mexican Hat etc. The choice of the filter is based onthe correlation between the shape of the signal with that of the filterresponse. The current shape of the PPG and ABP signal had a maximumcorrelation with the Daubechies filter with 4 vanishing moments (‘db4’).The flexibility of the design allows in choosing any filter and trainingthe neural networks based on that filter. In the current embodiment, achoice of db4 has been preferred.

The coefficients of the ‘db4’ or any chosen filter have high passcoefficients denoted by h_(l) and low pass coefficients denoted byg_(l). The signal is iteratively passed through the set of these filtersfor n times, where n ∈ {1,2 ...}. The total number of times, n, isdepended on the length of the signal in consideration. The presentembodiment considers a signal length of L sampled at F_(s) for both PPGand ABP signals. The number of decomposition level, n in present casecomes out to be:

n = ⌊log₂L⌋

The iterative decomposition of ABP and PPG using the circular filteroperation gives detailed high-pass components represented in equations(4) and (5).

$\psi\left( {PPG} \right)_{n,t} = {\sum\limits_{l = 0}^{L - 1}{{\widetilde{h}}_{nl} \times PPG_{t - l} \times mod(L)}}$

$\psi\left( {ABP} \right)_{n,t} = {\sum\limits_{l = 0}^{L - 1}{{\widetilde{h}}_{nl} \times ABP_{t - l} \times mod(L)}}$

Here, n is decomposition level, with n ∈ 1,2 ..., log2 L

-   L being the signal length-   h̃_(nl) being the coefficient of high pass filter periodised to    signal length L.

Similarly, the iterative approach yields the approximation componentsrepresented in equation (6) and (7)

$\phi\left( {PPG} \right)_{n,t} = {\sum\limits_{l = 0}^{L - 1}{{\widetilde{g}}_{nl} \times PPG_{t - l} \times mod(L)}}$

$\phi\left( {ABP} \right)_{n,t} = {\sum\limits_{l = 0}^{L - 1}{{\widetilde{g}}_{nl} \times ABP_{t - l} \times mod(L)}}$

g̃_(jl) being the coefficient of the low pass filter periodised to signallength L.

The relationship of the coefficients h̃_(nl) and g̃_(nl) with that of theselected filter which in this innovation was chosen to be ‘db4’ at everyiteration is modified through the scaling given by:

$\begin{array}{l}{{\widetilde{g}}_{nl} = \frac{g_{l}}{\sqrt{2^{n}}}} \\{{\widetilde{h}}_{nl} = \frac{h_{l}}{\sqrt{2^{n}}}}\end{array}$

The iterative decomposition provides vectors of Ψ(ABP)_(n,t) andΨ(PPG)_(n,t) which are significant for innovation. The low passcomponents of the MODWT decomposition, i.e. ϕ(PPG)_(n,t) andϕ(ABP)_(n,t) are ignored at each level are ignored. But in the n^(th)decomposition level of PPG and ABP, ϕ(PPG)_(n=n,t) and ϕ(ABP)_(n=n,t) isobtained, which is set to 0. This involves the low pass componentcontaining the mean of the signal. While PPG is demeaned resulting in 0energy at DC level, the ABP signal has a non-zero mean containing theMean Arterial Pressure whose estimation is described below.

The filter responses are determined for a signal length L chosen to be250. The responses demonstrate the changing bandpass filtering cut-offfrequency with every subsequent MODWT iteration with filter 1 having thehighest cut-off frequency at ƒ_(c) = ƒ_(s)/2² , filter 2 at fs/2³ and soon. The final iteration which is at n = log₂ 250 ≈ 7 yields a final lowpass component which is ignored.

The circular filter operation in equations (4) and (5) generates aconvolution between the periodised high pass filter coefficients at eachlevel with that of the signal giving a one-dimensional vector ofdecomposed signal MODWT coefficients. The periodisation converts thefilter coefficient at level n to a L × L square matrix giving the formof the matrix, where X represents either the PPG or ABP.

For n levels of decomposition, a 2-dimensional matrix for ABP and PPGare formed containing the high-frequency components,

$\psi(X) = \begin{bmatrix}{\psi_{1}\left( {X\lbrack t\rbrack} \right)} \\{\psi_{2}\left( {X\lbrack t\rbrack} \right)} \\ \vdots \\{\psi_{n}\left( {X\lbrack t\rbrack} \right)}\end{bmatrix}$

Data Structuring

The structuring of the data involved first the selection of thedecomposed coefficients is shown in FIG. 2 . The selection of thecoefficients was based on the energy level present at each level. Theidea is to incorporate those wavelet coefficients bands for PPG and ABP,which adds to the significant amount of energy. Since the estimationrequirement is that of ABP, the energy of the decomposed detailedcoefficients for ABP were analysed. Through the data analysis, it wasobtained to be in a signal length of L, the maximum density of energy isconcentrated in j number of bands, where j ⊂ n. Only those wavelet bandsare selected where the concentration of the energy adds up to 95% ormore.

In a wavelet decomposition of the signal length of L = 250 and n = 7,the values of j are obtained to be in levels 4-7. The bar graph in FIG.4 is evident of the levels 4-7 being the coefficients containing thehighest concentration of energy. Wavelets in these levels weresignificant while other remaining levels were zeroed, resulting in:

$\psi(X) = \begin{bmatrix}0 \\0 \\0 \\{\psi_{4}\left( {X\lbrack t\rbrack} \right)} \\{\psi_{5}\left( {X\lbrack t\rbrack} \right)} \\{\psi_{6}\left( {X\lbrack t\rbrack} \right)} \\{\psi_{7}\left( {X\lbrack t\rbrack} \right)}\end{bmatrix}$

Each resulting decomposed wavelet detailed coefficients of PPG,Ψ_(j)(PPG[t]) is organised into an autoregressive 2-D matrixrepresentation, Ψ_(j)PPG[t]). The dimension of the 2D matrix is ƒ × L.Here, ƒ ∈ Z⁺, is used as a set of features for the neural networks toconsider to estimate the present ABP Ψ_(n)[t = t_(present)] based onsome of the past and future samples. The generalised representation ofthe 2D matrix representation for a level n gives:

$\begin{array}{l}{\text{Ψ}_{j}\left( {PPG\lbrack t\rbrack} \right) =} \\\left\lbrack \begin{array}{lllll}{\psi_{j}\left( {PPG\left\lbrack {t - \frac{f}{2}} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t - \frac{f}{2} + 1} \right\rbrack} \right)} & \cdots & {\psi_{j}\left( {PPG\left\lbrack {t + L - \frac{f}{2} + 1} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + L - \frac{f}{2}} \right\rbrack} \right)} \\{\psi_{j}\left( {PPG\left\lbrack {t - \frac{f}{2} + 1} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t - \frac{f}{2} + 2} \right\rbrack} \right)} & \cdots & {\psi_{j}\left( {PPG\left\lbrack {t + L - \frac{f}{2}} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + L - \frac{f}{2} + 1} \right\rbrack} \right)} \\ \vdots & \vdots & \ddots & \vdots & \vdots \\{\psi_{j}\left( {PPG\lbrack t\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + 1} \right\rbrack} \right)} & \cdots & {\psi_{j}\left( {PPG\left\lbrack {t + L - 1} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + L} \right\rbrack} \right)} \\ \vdots & \vdots & \ddots & \vdots & \vdots \\{\psi_{j}\left( {PPG\left\lbrack {t + \frac{f}{2} - 2} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + \frac{f}{2} - 1} \right\rbrack} \right)} & \cdots & {\psi_{j}\left( {PPG\left\lbrack {t + L + \frac{f}{2} - 3} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + L + \frac{f}{2} - 2} \right\rbrack} \right)} \\{\psi_{j}\left( {PPG\left\lbrack {t + \frac{f}{2} - 1} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + \frac{f}{2}} \right\rbrack} \right)} & \cdots & {\psi_{j}\left( {PPG\left\lbrack {t + L + \frac{f}{2} - 2} \right\rbrack} \right)} & {\psi_{j}\left( {PPG\left\lbrack {t + L + \frac{f}{2} - 1} \right\rbrack} \right)}\end{array} \right\rbrack\end{array}$

From each signal of length L, wavelet decomposed to levels n andconsidering the highest energy concentration of the wavelet coefficientsonly in levels j, gave length of j equivalent number of 2-D matrices.

The output ABP for each significant level, j was arranged in vectorizedform giving,

$\begin{array}{l}{\psi_{j}\left( {ABP\lbrack t\rbrack} \right) =} \\\left\lbrack {\psi_{j}\left( {ABP\lbrack t\rbrack} \right),\psi_{j}\left( {ABP\left\lbrack {t + 1} \right\rbrack} \right),\mspace{6mu}\ldots,\psi_{j}\left( {ABP\left\lbrack {t + L} \right\rbrack} \right)} \right\rbrack\end{array}$

LSTM Neural Network

Long short-term memory (LSTM) is a variant of neural networks. They haveinternal gates. These gates ensure the information it is obtaining fromthe previous values are selectively read, selectively forgotten andselectively written. This selective process makes the LSTM efficient inremembering only crucial information from the past rather thaneverything.

The structure of a neuron of an LSTM consist of forgetting, read, writeand state gates. These gates have associated parameters - weights andbiases, which are learnt through the process of feedforward andback-propagation. The size of the weights and biases are depended on thedesigner of the neural networks and the purpose for which it isdesigned. In the present innovation, the single LSTM neuron has beenassumed to have α units. These units are the matrix sizes of theparameters inside which have to be learnt in the training process. Thereare 12 essential parameters associated with the 4 essential gatesrequired to be learnt. These include the recurrent weights given byW_(x), kernel input weights are given by U_(x) and biases b_(x).

$Parameters = \begin{pmatrix}{W_{i},U_{i},b_{i}} \\{W_{f},U_{f},b_{f}} \\{W,U,b} \\{W_{o},U_{o},b_{o}}\end{pmatrix}$

The neuron structure of an LSTM is shown in FIG. 5 . The neuron receivesinformation from the previous time-step in the form of the estimatedoutput h_(t-1) and the cell state C_(t)-₁. The Previous Cell Statecontains the information of factors such as what has been remembered,read and written in the previous state. It may also include what hasbeen forgotten in the previous state. The Past Output h_(t-1) is theoutput from the previous time-step. In its simplistic form, it is thepredicted output of the wavelet coefficient of Ψ_(j)(4BP[t - 1]).

The input to a single neuron as shown in FIG. 5 is the column vector attime t of 2D matrix Ψ_(j)(PPG[t]) of level j. The column vector contains

$\frac{f}{2}$

past samples and

$\frac{f}{2} - 1$

future samples from present time state t. The vector input along withgates estimate the output h_(t).

Internally the first gate is the forget gate. The information it isacquiring from the previous state is modified to incorporate theinformation from the past and forget unwanted or less involvedinformation. The output range of f_(t) is between 0 and 1.

The next gate is the input which decides what information from thepresent input Ψj(PPG[t]) must be read. This selective reading isacquired through the process of learning to enable it to have selectivefeature selection extraction.

The third cell stage is the temporary cell state, c̃_(t). The temporarycell state is used for the modification of the previous cell state,incorporating the information from the present input as well.

Finally, the output, o_(t) ensures selective writing in predicting thepresent output.

The present cell state, C_(t) and present cell output h_(t) is given by:

$\begin{matrix}{C_{t} = f_{t} \times C_{t - 1} + i_{t} \times {\widetilde{c}}_{t}} \\{h_{t} = o_{t} \times \tanh C_{t}}\end{matrix}$

A recursive dependency is observed from (14), with the output dependingon the different gates, which in turn depends on the previous states.

Neural Network Structure and ABP Coefficients Estimation

The overall design of the expected neural network is shown in FIG. 6 .This is a representation of the deep LSTM which is used for theestimation of a single wavelet coefficient at level j from thestructured matrix of PPG wavelet at that level. A multi-layer LSTM ispreferred to estimate the non-linear patterns in the dynamic nature ofthe body.

Each layer of the LSTM neural network follows L number of cells based onthe signal length L. This spans over a depth of i consisting ofdifferent structural layers. The flexibility of the design allowschanging the value of i based on the accuracy of the system. But largervalues of i would require higher computation power and may causeoverfitting leading to inaccuracy in the actual result. Currentinnovation involves the use of i to be 3 layers.

In the current innovation, at each layer, the number of units presentinside each neuron is varied which changes the number of Parameters tobe learnt to develop a complex time-based correlation. Further, theinitialisation of the LSTM states C_(t) and h_(t) initially at everycell level was set to empty matrix 0.

The final deep LSTM layers are followed by a single deep layer ofartificial neural network. This is used to combine the estimates of theoutput matrices h_(t) in the LSTM layer i and provide a single waveletpressure estimate Ψj(ABP[t]). The subsequent ABP wavelet estimates are acomplex combination of the past estimates in every layer giving a finaloutput of a vector consisting of Ψ_(j)(ABP[t]).

In the training phase, a standard back-propagation technique is used inlearning the Parameters. This is done through the feedback of the outputerror generated at every epoch on the training data (obtained from theMIMIC-II database) with the ground truth pressure wavelets andpropagating the error backwards. This backpropagation of the errorresults in the update of the weights, as shown in FIG. 3 . Furthermore,based on the number of desired levels under consideration possessing asignificant amount of energy, that many numbers of LSTM neural networkswere trained, which in this case is the number of levels, j. In apreferred implementation, the number of structures trained was 4 with asignal length of L = 250.

Once the weights of all the layers are learnt, the implementation phaseinvolves structuring and presenting the inputs to the network, whichaccumulates and multiplies the values to obtain the final estimatedoutput vector.

Reconstruction of Pulse Pressure Waveform

The reconstruction of the pulse pressure waveform follows the outputvectors generated by the length of j number of LSTM for each level. Theoutput vector contains

$\psi\left( {ABP\lbrack t\rbrack} \right) = \begin{bmatrix}{\psi_{j{(1)}}\left( \overline{ABP\lbrack t\rbrack} \right)} \\{\psi_{j{(2)}}\left( \overline{ABP\lbrack t\rbrack} \right)} \\ \vdots \\{\psi_{j{({len{(j)}})}}\left( \overline{ABP\lbrack t\rbrack} \right)}\end{bmatrix}$

Knowing the fact that the reconstruction is depended on the levels whichhave been considered in the LSTM models were significant, and the levelswhich were not included in the LSTM estimation especially thecoefficients lying in higher frequency bandwidth were zeroed, forming a2D vector of size (log₂ L) × L,

$\psi\left( {ABP\lbrack t\rbrack} \right) = \begin{bmatrix}0 \\0 \\0 \\{\psi_{j_{1}}\left( {X\lbrack t\rbrack} \right)} \\{\psi_{j_{2}}\left( {X\lbrack t\rbrack} \right)} \\{\psi_{j_{3}}\left( {X\lbrack t\rbrack} \right)} \\{\psi_{j_{4}}\left( {X\lbrack t\rbrack} \right)}\end{bmatrix}$

The reconstruction filter involves the use of complementary filters ofthe analysis decomposition filter. In the current innovation, usingDaubechies filter, ‘db4’ and applicable to any other wavelet filters inthe family of wavelets, the relationship between the h_(l) which is thedecomposition filter with that of hr_(l).

The inverse MODWT reconstruction is obtained by:

$X(t) = {\sum\limits_{j = 1}^{n}{\sum\limits_{l = 0}^{L - 1}{\widetilde{hr}j,_{l}\psi_{j}\left( {X\lbrack t\rbrack} \right)}}} + {\sum\limits_{l = 0}^{L - 1}{{\widetilde{gr}}_{j,l}\phi_{j = j_{o}}\left( {X\lbrack t\rbrack} \right)}}$

In equation (17), the latter part involving the low pass component withϕj is insignificant in the present embodiment as the DC components havenot been involved in the reconstruction using LSTM. Hence equation (17)reduces to:

$\begin{matrix}{X(t) = {\sum\limits_{j = 1}^{n}{\sum\limits_{l = 0}^{L - 1}{{\widetilde{hr}}_{j,l}\psi_{j}\left( {X\lbrack t\rbrack} \right)}}}} \\{X(t) =} \\{{\sum\limits_{l = 0}^{L - 1}{{\widetilde{hr}}_{j_{1},l}\psi_{j_{1}}\left( {X\lbrack t\rbrack} \right)}} + {\sum\limits_{l = 0}^{L - 1}{{\widetilde{hr}}_{j_{2},l}\psi_{j_{2}}\left( {X\lbrack t\rbrack} \right)}} + {\sum\limits_{l = 0}^{L - 1}{{\widetilde{hr}}_{j_{3},l}\psi_{j_{3}}\left( {X\lbrack t\rbrack} \right)}}} \\{+ {\sum\limits_{l = 0}^{L - 1}{{\widetilde{hr}}_{j_{4},l}\psi_{j_{4}}\left( {X\lbrack t\rbrack} \right)}}}\end{matrix}$

Following equation (18), in the expanded form, at every level, thefilter coefficients

${\widetilde{hr}}_{j,\, l}$

is modified following the representation in equation (8),

${\widetilde{hr}}_{j,\, l} = \frac{hr_{l}}{\sqrt{2^{j}}}$

Leading from the expanded form in equation (18), the MODWT calculationleads to obtaining the pulse pressure, ΔP = X(t).

One embodiment of a method 100 of estimating blood pressure of a subjectis summarised in FIG. 7 . The method 100 shown in FIG. 7 includes thesteps of: receiving 102 a photoplethysmogram (PPG) signal from alight-based Pulse-Plethysmography sensor applied to the skin of thesubject; processing 104 the PPG signal using a Wavelet transformalgorithm to derive PPG wavelet coefficients in a plurality of frequencybands; estimating 106 blood pressure coefficients by processing the PPGwavelet coefficients using a machine learning algorithm that has beentrained on training data including PPG wavelet coefficients derived fromPPGs from light-based Pulse-Plethysmography sensors applied to the skinof test subjects correlated with Invasive Arterial Blood Pressure (ABP)wavelet coefficients derived from simultaneously obtained ABPmeasurements of systolic and diastolic blood pressure of the testsubjects; and reconstructing 108 a pulse blood pressure waveform betweensystolic and diastolic blood pressure of the subject from the bloodpressure coefficients.

Another embodiment of the method 100 is shown in FIG. 8 in respect ofthe further method 110 of estimating blood pressure of a subject. Thefurther method 110 includes the steps of: processing 112 the pulse bloodpressure waveform using a further Wavelet transform algorithm to derivepulse blood pressure wavelet coefficients in a plurality of frequencybands; extracting 114 features from the pulse blood pressure waveletcoefficients in the plurality of frequency bands; estimating 116 meanarterial blood pressure (MAP), systolic blood pressure (SBP) ordiastolic blood pressure (DBP) coefficients by processing the featuresusing a further machine learning algorithm that has been trained on thetraining data; reconstructing 118 the MAP, SBP and or DBP waveforms fromthe MAP, SBP and or DBP coefficients, respectively; and combining 120the pulse blood pressure waveform and one or more of the MAP, SBP or DBPwaveforms to reconstruct an absolute blood pressure waveform of thesubject.

The absolute blood pressure waveform of a human subject is shown in FIG.9 . The maximum and minimum pressure are the Systolic (SBP) andDiastolic (DBP) blood pressures, respectively, which forms the PulsePressure Waveform. The mean of this waveform is the Mean ArterialPressure (MAP). Embodiments of the present invention estimate the SBP,DBP and MAP waveforms over time using deep learning models. It cantherefore be seen that, by combining the MAP, SBP or DBP waveforms,having pressure values, with the pulse pressure waveform, an absoluteblood pressure waveform of the subject can be derived.

In respect of this method 110, the hardware used for the collection ofbio-signals from a subject is a modified version of aPulse-Plethysmography (PPG) sensor operating at a wavelength of 805 nm,also called the isosbestic wavelength, as above.

In an embodiment of this method 110 of estimating blood pressure, theacquired waveform from the PPG sensor is passed through a set ofalgorithms comprising: a Maximally Overlap Discrete Wavelet Transformalgorithm; an LSTM Based Pulse Pressure Generation Algorithm; a FeatureMatrix Generation algorithm; and a Conv-LSTM and Deep Neural NetworkBased Pressure Estimation algorithm.

FIG. 10 shows a summary of an embodiment of the method 110 of estimatingblood pressure waveform from a combination of the abovementioned LSTMderived pulse pressure waveform algorithm and MAP waveform estimationalgorithm.

In an embodiment, the method 100 and the further method 110 use hardwareshown in FIG. 11 . The hardware includes multiple LEDs and Photodiodesfor acquisition of the PPG signal from any part of a subject’s wrist.The hardware further includes a micro-controller to control the sensorswith Bluetooth and WIFI capabilities for data transfer; apulse-plethysmography PPG sensor chip. The multiple LEDs andPhoto-diodes operate at wavelength range between 800 - 850 nm(isosbestic) and the hardware further includes multiplexers controllingthe different pairs. In addition, the hardware includes an accelerometerto remove any motion artefacts, power supply units and a batterycharging controller.

Multiple pairs or groups of LEDs and Photodiodes ensures that datacollection is from any location of the wrist or limb of the subjectwhere the device is placed. FIG. 12 shows a block diagram of thehardware in use.

A master micro-controller controls the pair of LED and photodiode toturn on and off based through an N-input multiplexer. The activated pairof the LED-photodiode is controlled by the PPG Sensor RX-TX controlchip. This chip drives the LED at specific pulse-repetition frequency,drive voltage and current. It receives a signal from the selectedphotodiode. The signal is filtered and amplified. The processed signalis sent into the micro-controller which then decides the SNR of thesignal.

The SNR of the signal will be affected by the motion artefact present orthe location at which the pair of LED and photodiode is present. Themotion artefact is calculated using the on-board 6 Axis IMU. Themicrocontroller runs the algorithm to select the most appropriate pairof the LED and photodiode, as embodied in the flow chart of FIG. 12 .That is, the microcontroller calculates the SNR of each LED-photodiodepair along with a compensation value from the IMU sensor(s). Themicro-controller than selects the pair of LED-photodiode with the bestSNR.

An embodiment of the method of estimating blood pressure in use isdescribed with reference to the flow chart of FIG. 13 . As mentioned,the PPG sensors are formed of multiple pairs of LEDs and Photodiodes,which are operating at the isosbestic wavelength. The sensors arepreferably disposed in a device that sits around the wrist and thesensors are equally spaced. This ensures that any motion of the deviceworn of the wrist does not hinder data collection. Further, it willensure that the best Signal-To-Noise (SNR) PPG signal is collected andused, which will be highly dependent on the location of the radial orbrachial artery on the subject and where the device is intended to besitting on the subject.

In an initial cycle, all the pairs of LEDs and Photodiodes will beactivated one after the other. A small batch of signals of a few secondswill be collected from each pair. The collected signal from each pair isprocessed to see the presence of:

-   sinusoidal oscillations matching the centre frequency equivalent to    heart rate;-   power spectrum analysis and comparison - higher powers around the    desired frequency regions for better SNR; and-   lower power for harmonics, if the light waves are reflected from    bones or tissues modifying or interfering with the centre frequency.

The pair of the LED and Photodiode adhering to the selection criteriawill be considered to have the best SNR. The signal from this pair willbe used subsequently for purpose of the estimation.

FIG. 13 further shows the signal from the PPG sensor being transmittedusing a wireless method to a device which can undertake machine learningestimation, processing and reconstruction. The transfer can be throughBluetooth or WIFI. This is flexible and the mode of communication can bechanged based on the target device, such as a mobile phone operating adedicated application or a custom designed monitor.

FIG. 13 further includes the step of Pre-processing the signal from thePPG sensor. The initial pre-processing is followed from thereconstructing the pulse blood pressure waveform, which involves lowpass filtering with a cut-off frequency at 20 Hz. In this frequencyrange, the expectation of higher physiological response is present.

In the LSTM-MODWT Based 4P step, the pre-processed signal is passedthrough the LSTM based 4P estimation algorithm which has been describedabove.

In the Feature Extraction step, the feature generation process involvesautomatic extraction of 2D features estimated from the LSTM algorithm,hand-designed specific features from PPG signal and hand-designedspecific features from the estimated 4P waveform. The three sets offeatures are processed, combined and concatenated to form the trainingand prediction algorithm for the method of estimating absolute bloodpressure.

In the ConvLSTM and DNN algorithms step, the three sets of features fromthe feature extraction process are processed. The automatic 2D featuresare passed into the ConvLSTM algorithm while the hand-designed featuresare passed into multi-layer perceptrons (MLP). This results in formingthree branches for the network getting three different sets of inputs.Subsequently, the algorithm ensures the concatenation of the threebranches to form a convoluted semi-supervised structure. The weights ofthe structure is learnt using the standard Back-Propagation algorithm.The final output of the structure has been designed to be either theSBP, DBP and MAP depending on what is being learnt. The optimisation ofthe structure with different hyper-parameters is based on the finaloutput. The appropriate hyper-parameters are selected based on theperformance for these three sets of the outputs.

In the Rescaling and Median Filtering step, following the estimationfrom the previous step, the output, be it SBP, DBP or MAP, is set todown-sampled vectors consisting N number of samples. In the step, N hasbeen taken to be ⅒th the length of the actual input signal. A medianfiltering is carried out for every subsequent N length of the outputSBP, DBP or MAP ensuring that any sudden changes in the estimatedpressure output due to errors in the networks are minimised.

In the final Absolute Total Pressure Waveform step, the estimatedwaveforms of MAP, SBP and DBP are added to, or combined with, the ΔPwaveform estimated and reconstructed according to the above described togenerate the absolute (and continuous) blood pressure waveform of thesubject.

The method 110 is described in more detail below. The method 110commences with a training phase. As above, Invasive Arterial BloodPressure (ABP) measurements and PPG signals are recorded from manypatients simultaneously. As described with respect to FIG. 3 , the PPGand ABP signals are passed through the Maximally Overlap DiscreteWavelet Transform (MODWT) algorithm and broken down into differentbandwidths of frequencies between 0.5-8 Hz. The PPG MODWT waveformbetween specific frequency bandwidths is created into auto-regressivestructures which are correlated to the ABP MODWT waveform between thatspecific bandwidth. The error is calculated after every iteration oftraining and the LSTM weights are updated.

The final result from this is the complex correlation time seriesnetworks relating the MODWT coefficients of PPG and ABP at differentbandwidths as shown in the example FIGS. 19 a-19 c . The estimatedwaveforms shown in these FIGS. 19 a-19 c as subplot 2,3,4 and 5represent the wavelet detailed coefficients which were estimated usingthe LSTM-MODWT at bandwidth 4-8 Hz (ψ4), 2-4 Hz (ψ5), 1-2 Hz (ψ6) and0.5-1 Hz (ψ7). These outputs from the different bandwidths are combinedusing the inverse MODWT to reconstruct the pulse pressure waveform, ΔP.

The reconstructed pulse pressure waveform, ΔP, along with the PPG signalare used as input features for the method 110. The overall design of thetraining model of the method 110 is shown in FIG. 14 . In particular,FIG. 14 shows the training phase of a Conv-LST-DNN based AbsolutePressure Estimation Model.

The estimated bands of the detailed coefficients are arranged intotime-frequency 2D matrices which are used as features for the ConvLSTMnetwork. The inverse MODWT reconstructed ΔP and PPG waveform are used toderive some of the hand designed features involving the time andamplitude values. These are passed into two separate branches ofmulti-layer perceptrons. All the branches are concatenated and the finaloutput the expected estimate of the absolute pressure which may beeither SBP, DBP or MAP. The error is calculated with the ground truthpressure values and the error is backpropagated to update the neuralnetwork weights.

Convolutional-LSTM Network

Convolutional-LSTM (ConvLSTM) is a hybrid of Convolutional NeuralNetworks (CNN) and Long Short Term Memory (LSTM) Recurrent NeuralNetworks. LSTMs have the capability of temporal forecasting. They aredesigned to remember the long and short state of the system, making thembeneficial for time-series forecasting.

ConvLSTM consists of cells as shown in FIG. 15 . The structure of aneuron of a ConvLSTM consists of forgetting, read, write and state gatesbut, instead of these being 1D vectors, they are 2D matrices. Thesegates have associated parameters - weights and biases - but a depthparameter is added as well with the number of filters similar to a CNNwhich are learnt through the process of feedforward andback-propagation. The size of the weights and biases are dependent onthe designer of the neural networks and the purpose for which it isdesigned. Preferably, the single ConvLSTM neuron has been assumed tohave α units which are kept variable based on the depth and performance.These units are the matrix sizes of the parameters inside each cellwhich have to be learnt in the training process.

In the 2D form, the equations for the gates and the output is shown inequations (20)-(25), and the neuron shown in FIG. 15 . All themultiplications have been replaced with the convolution between theinputs and the multiple number of filters.

f_(t) = σ(W_(f) ⊗ h_(t − 1) + U_(f) ⊗ X_(t) + b_(f))

i_(t) = σ(W_(i) ⊗ h_(t − 1) + U_(i) ⊗ X_(t) + b_(i))

c_(t) = tanh(W ⊗ h_(t − 1) + U ⊗ X_(t) + b)

o_(t) = σ(W_(o) ⊗ h_(t − 1) + U_(o) ⊗ X_(t) + b_(o))

$C_{t} = C_{t - 1} \otimes f_{t} + i_{t} \otimes \widetilde{c_{t}}$

h_(t) = o_(t) ⊗ tanh C_(t)

where,

-   all the variables are 2D,-   f_(t) = forget gate-   i_(t) = selective read gate-   C_(t) = state gate-   o_(t) = selective write gate-   W= corresponding recurrent weights of the gates-   b = corresponding biases of the gates-   U = corresponding input weights of the gates-   X_(t) = 2D input matrix-   C_(t-1) = previous cell state-   h_(t-1) = previous cell-   h_(t) = current cell output

The ConvLSTM neuron as shown in FIG. 15 takes a 2D matrix of input atspecific time index and passes it through four gates. But unlike LSTMwhich consists of four 1D weight vectors each of length n, the ConvLSTMconsists of n̅ 2D filters, which makes it 4 × n̅ 2D filters to be learntfor each neuron during the training process.

Passing of 1 feature matrix, (X_(t)) of dimension D₁ × H₁ through a n̅number of 2D filters in a gate of shape K × K ∈ {K = N} results in 3Dtensors of shape of n̅ × D₂ × H₂. These are the number of featuresextracted for each gate after convolution, at a timeframe of t = t_(n).This has a capacity to derive extremely deep and abstract features whichwould have been missed in a simple neural network. The shape of D₂ andH₂ highly depends on the hyperparameters selected for the ConvLSTM,especially the dimension of filter K, stride of the filter, zero paddingand number of filters n̅. Usually, keeping the filter dimension K = 3,stride as 1, zero padding as 1 (which refers to padding the boundariesin order to preserve the same shape of the output as that of the input)and number of filters as 64, the output dimension from each gate wouldcalculate to be:

$D_{2} = \frac{D_{1} - K + 2 \times Padding}{Stride} + 1 = 4$

$H_{2} = \frac{H_{1} - K + 2 \times Padding}{Stride} + 1 = \left( {H_{1} - 3 + 2} \right) + 1 = H_{1}$

Hence, with this example one gate would result in 64 × 4 × H₁ number offeatures. Here, H₁ depends on the length of the feature vector in timedomain. Undertaking the operations from each gates having the dimensionof 64 × 4 × H₁ shown in equations (20) - (25) would result in h_(t) andC_(t) as the output from the neuron of dimensions 64 × 4 × H₁.

In a multi-layer structure, the number of time inputs depends on thenumber of frames the inputs have been divided. Considering each frame oflength X, then X number of ConvLSTM cells will be stacked sequentially,similar to the LSTM network as shown in FIG. 16 . This further wouldlead to X number of output matrices for h_(t).

The output from each neuron has features estimated from the 2D matrix,which when passed through the batch normalisation and subsequentlypooling layer are reduced in spatial shape. In the final step, theoutput from X frames are flattened and passed through a dense layer as avector to run the regression analysis. The dense layer is called a TimeDistributed layer. This takes the flattened inputs from each time stepand generate a vector of time sequence output as shown in the FIG. 17 .

Data Structuring and Features

The output from the n LSTM networks are n vectors of length L. In anembodiment, n = 4 and L = 250. n represents the number of bandwidths theLSTM outputs which are divided into 4 regions consisting of: 0.5 - 1 Hz,1 - 2 Hz, 2 - 4 Hz and 4 - 8 Hz. The length of the signal chosen as thisis based on the sampling frequency, which in this case is F_(s) = 125Hz. This captures 2 second of the PPG signal and the corresponding ΔPfrom the LSTM output. The n number of vectors are concatenated into n ×L matrix,

$\psi^{ABP}\lbrack k\rbrack = \begin{bmatrix}0 \\0 \\0 \\{\psi_{1}\left\lbrack {X\lbrack t\rbrack} \right\rbrack} \\{\psi_{2}\left\lbrack {X\lbrack t\rbrack} \right\rbrack} \\{\psi_{3}\left\lbrack {X\lbrack t\rbrack} \right\rbrack} \\{\psi_{4}\left\lbrack {X\lbrack t\rbrack} \right\rbrack}\end{bmatrix}$

where:

-   ψ(ABP[t]) = log₂L × L matrix of the decomposition wavelet    coefficients-   ψ_(n)(.) = vector output from the n^(th) LSTM models-   X(t) = feature matrix for LSTM at particular bandwidth of frequency

The presence of the zero padding are due to the higher frequencydomains, greater than 16 Hz, which have not been considered to be regionof interest of the PPG-blood pressure correlation.

In order to reconstruct the signal, the inverse MODWT of the signal inequation (28) was taken. As discussed above, the output was passedthrough a Daubechies filter with 4 vanishing moments, db4. Eachestimated vector is convolved with the filter and added in order toreconstruct the pulse pressure. A simplified version of the equationprocess is shown in the equation (29).

$\text{Δ}P(t) = {\sum\limits_{i = 0}^{n}\left\lbrack {\text{ψ}\left( {ABP\lbrack t\rbrack} \right)_{t} \otimes \widetilde{h_{r}}} \right\rbrack}$

where:

-   ΔP(t) = pulse pressure-   h̃_(r) = Daubechies reconstruction filter coefficients

In the method 110, the ψ(ABP [t]) equation (28) as well as the ΔP(t)equation (29) were used in the generation of the feature matrix.

ConvLSTM cells accept 2D matrices as inputs. In this case, the estimatednon-zero MODWT coefficients of ABP as shown in equation (28) werestructured into inputs. The outputs on the other hand were theapproximation coefficients of the ABP expected to be estimated from themachine learning algorithm.

The estimated detailed wavelet coefficients following equation (28) ofsignal of length L comes out to be of shape 4 × N shown in equation(30). Where, as discussed, the estimated sequence-to-sequencerelationship was developed for four levels consisting of highest energyfrom above.

$\begin{array}{l}{\left\lbrack \begin{array}{l}{\text{ψ}_{1}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)} \\{\text{ψ}_{2}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)} \\{\text{ψ}_{3}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)} \\{\text{ψ}_{4}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)}\end{array} \right\rbrack =} \\\left\lbrack \begin{array}{llll}{\text{ψ}_{1}\left( {ABP\lbrack 0\rbrack} \right)} & {\text{ψ}_{1}\left( {ABP\lbrack 1\rbrack} \right)} & \cdots & {\text{ψ}_{1}\left( {ABP\lbrack L\rbrack} \right)} \\{\psi_{2}\left( {ABP\lbrack 0\rbrack} \right)} & {\text{ψ}_{2}\left( {ABP\lbrack 1\rbrack} \right)} & \cdots & {\text{ψ}_{2}\left( {ABP\lbrack L\rbrack} \right)} \\{\psi_{3}\left( {ABP\lbrack 0\rbrack} \right)} & {\psi_{3}\left( {ABP\lbrack 1\rbrack} \right)} & \cdots & {\psi_{3}\left( {ABP\lbrack L\rbrack} \right)} \\{\psi_{4}\left( {ABP\lbrack 0\rbrack} \right)} & {\psi_{4}\left( {ABP\lbrack 1\rbrack} \right)} & \cdots & {\psi_{4}\left( {ABP\lbrack L\rbrack} \right)}\end{array} \right\rbrack\end{array}$

The 2D matrix in equation (30) was structured into 4 × M dimensionalsmaller 2D feature matrices with non-overlapping parts. Here M is thefraction of signal length L. This further resulted in L/M number oftotal feature matrices for a signal. Mathematically one feature matrixis shown in equation (31).

$\begin{array}{l}{\left\lbrack \begin{array}{l}{\psi_{1}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)} \\{\psi_{2}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)} \\{\psi_{3}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)} \\{\psi_{4}\left( \overline{\overline{ABP\lbrack t\rbrack}} \right)}\end{array} \right\rbrack =} \\\left\lbrack \begin{array}{llll}{\psi_{1}\left( {ABP\lbrack t\rbrack} \right)} & {\psi_{1}\left( {ABP\left\lbrack {t + 1} \right\rbrack} \right)} & \cdots & {\psi_{1}\left( {ABP\left\lbrack {t + M} \right\rbrack} \right)} \\{\psi_{2}\left( {ABP\lbrack t\rbrack} \right)} & {\psi_{2}\left( {ABP\left\lbrack {t + 1} \right\rbrack} \right)} & \cdots & {\psi_{2}\left( {ABP\left\lbrack {t + M} \right\rbrack} \right)} \\{\psi_{3}\left( {ABP\lbrack t\rbrack} \right)} & {\psi_{3}\left( {ABP\left\lbrack {t + 1} \right\rbrack} \right)} & \cdots & {\psi_{3}\left( {ABP\left\lbrack {t + M} \right\rbrack} \right)} \\{\psi_{4}\left( {ABP\lbrack t\rbrack} \right)} & {\psi_{4}\left( {ABP\left\lbrack {t + 1} \right\rbrack} \right)} & \cdots & {\psi_{4}\left( {ABP\left\lbrack {t + M} \right\rbrack} \right)}\end{array} \right\rbrack\end{array}$

Here, n ∈ 1, 2, ..., L/M. In the embodiment, M is taken to be 10. Thisresults in providing 10 samples of the wavelets approximationcoefficient matrix.

The output vector is the approximation coefficients of the signal oflength L or the bandwidth of the signal between 0-0.5 Hz represented byϕ. This ϕ can also be either the upper or the lower enveloperepresenting the SBP or DBP waveform of the blood pressure. This vectoris of length L. The approximation coefficients vector was sub-dividedinto L/M output vectors. With the input feature being of dimension 4 ×M, the output was selected to be the median location of 1 × Mapproximation coefficients. This was done in order to incorporate anon-causal model with the output depending on a few past samples and afew future samples shown in equation (32).

$y = \phi\left( {ABP_{\frac{t_{n} + M}{2}}} \right)$

where t_(n) is the starting location of each of the feature matrices

Hence, the output vector from the ConvLSTM is expected to be an N/Mlength of approximation coefficients shown in equation (33).

$y = \begin{bmatrix}{\phi\left( {ABP_{\frac{t_{1} + M}{2}}} \right)} & {\phi\left( {ABP_{\frac{t_{2} + M}{2}}} \right)} & \cdots & {\phi\left( {ABP_{\frac{t_{\frac{N}{M}} + M}{2}}} \right)}\end{bmatrix}$

A diagrammatic representation is shown in FIGS. 19 a-19 c . Theestimates ψ₄ - ψ₇ are the estimated wavelet bands as described above.These are concatenated in a 2D form. A window moves over creating acorrelation with the ϕ approximation coefficients vector at the medianpoint as shown in FIGS. 19 a-19 c . This entire process occurs withinthe ConvLSTM network.

Deep Neural Network (DNN) Features

In an embodiment, hand-designed features are selected for further bloodpressure estimation. These selected features selected were heuristicwith the number of features ranging from 5 to 50, all extractedmanually. The features include both common features based on time andamplitude, such as rise time, fall time, rising inflection time, fallinginflection time, area under the rising curve, area under the fallingcurve and heart rate (centre frequency of the signal), and features fromthe current approximation coefficient estimation algorithm. By usingthese features, the accuracy of the estimation of blood pressure isimproved.

But the addition of these features pose a problem to the ConvLSTMNetwork as they take time-series 2D matrices of the detailedcoefficients. One possibility was training of different networks, onebased on the ConvLSTM architecture and one based on the deep networkswith hand-designed features. This is the method of Ensemble Learning.The estimated accuracy will be the weighted average of the two differentmodels. The problem with this technique is the learned weights of eachmodel will be independent from each other. As the algorithm has to sharethe output error from the back propagation, a different architecture wasselected to be used.

In this embodiment, multiple models were used which are concatenated atthe end of the process to estimate a single output. This is based on thetheory of developing a functional architecture rather than a sequentialarchitecture. In a Sequential Architecture, the data flows from theinput to output without branching, while in a Functional Architecture,branches are present making it possible to further add different inputswhich have a possibility of improving the network accuracy. In thisembodiment, functional architecture is used.

A deep neural network was used in this case which has been designed toincorporate selected (hand-designed) features from:

-   PPG waveform; and-   Pulse pressure waveform.

Other important point to consider is the absolute values of time andamplitudes will vary significantly among the signals. Instead of takingthese absolute values, ratios were obtained with respect to a totallength of the window of the signal in cases of time related features andrange of the signal in case of amplitude related the features.

A total of 13 features from PPG were obtained which are listed below.These features were averaged for a signal with N heart cycles of PPG:

-   1. Total Window Time (WT) - time duration of one entire cycle of the    PPG signal.-   2. Rise Time Ratio (RTR) - the ratio of the time taken for the PPG    to rise to the maximum value to that of the WT.-   3. Fall Time Ratio (FTR) - the ratio of the time taken for the PPG    to fall from the maximum to minimum value to that of the WT.-   4. Rise Inflection Time Ratio (RITR) - the ratio of the time taken    to reach the first rising inflection point to that of the WT.-   5. Fall Inflection Time Ratio (FITR) - the ratio of the time taken    to fall from the maximum value to the falling inflection point to    that of the WT.-   6. Heart Rate Frequency (HRF) - the instantaneous heart rate    calculated from the N heart cycles.-   7. Pulse Plethysmogram Intensity Ratio (PIR) - the ratio of the    maximum to the minimum value of the PPG signal.-   $PIR\mspace{6mu} = \mspace{6mu}\frac{V_{PPG_{MAX}}}{V_{PPG_{MIN}}}$-   8. Rising Inflection Point Relative Height (RIPRH) - the ratio of    the PPG amplitude at the rising inflection point to the maximum    height.-   9. Falling Inflection Point Relative Height (FIPRH) - the ratio of    the PPG amplitude at the falling inflection point to the maximum    height.-   10. Systolic Area Ratio (SAR) - the area under the rising part of    the PPG signal to the total area under the curve.-   11. Diastolic Area Ratio (DAR) - the area under the falling part of    the PPG signal to the total area under the curve.-   12. Systolic to Diastolic Area Ratio (SDAR) - the ratio of the SAR    to DAR.-   13. AC to DC Ratio (ADR) - the ratio of the PPG amplitude to the RMS    of the signal:-   $ADR\mspace{6mu} = \mspace{6mu}\frac{PPG_{max} - \mspace{6mu} PPG_{min}}{rms\left( {PPG} \right)}$

From the Pulse Pressure, a total of 7 features related to the amplitudewere selected listed below:

-   1. Maximum of the pulse pressure-   2. Minimum of the pulse pressure-   3. Amplitude of the pulse pressure-   4. Systolic area ratio-   5. Diastolic area ratio-   6. Rising inflection height-   7. Falling inflection height.

Integration - ConvLSTM and DNN Networks

The ConvLSTM and the DNN networks have been designed in a functionalstructure with three branches. One branch taking the input as the timeseries wavelet matrices shown in equation (31). The second branch takesthe PPG features as the input while the third branch takes the pulsepressure features as input.

The output is the estimation of the approximation coefficient vectorshown in equation (33). The overall network structure with the branchedinputs are further shown in FIG. 17 .

The output from the time distributed dense of the ConvLSTM isconcatenated with the outputs of the two DNN in the functional networkstructure. The concatenated output is passed through multiple layers ofhidden neurons in Output Network. Finally, the loss is calculatedbetween the ground truth approximation coefficients and the predictedapproximation coefficients of the batch during training andbackpropagated with the error.

The final trained model is expected to estimate a vector sequence oflength L/M number of approximation coefficients ϕ. Flexibility of thedesign also allows to train the model to estimate the vector output forSBP or DBP instead of the ϕ which represents the MAP.

Post Processing and Median Filtering

The final output ϕ is a vector L/M, which in the current case has beentaken to be

$\left\lfloor \frac{250}{10} \right\rfloor\mspace{6mu} - \mspace{6mu} 1\mspace{6mu} = \, 24$

. The negative 1 is added in the current case to have an even samplelength which can be used in creating a moving median filter. Instead ofgiving a capture of MAP, SBP or DBP - based on what the ϕ is learnt tobe, a median filtering is carried out by estimating on the subsequentsignal of length L.

In an embodiment, L is taken to be of length of 2 seconds (samplingfrequency = 125 Hz). If a signal of length 10 seconds is collected fromthe PPG, the signal is broken down into 40 overlapping frames 2 secondsworth of signal, with 75% overlap. This overlap can be more or less.Each sequence is used for the estimation of ϕ vector, resulting in 40vectors in the current instance.

A moving median filter is applied across the overlapped signal. Themedian filter has the advantage of being less prone to outliers ascompared to the moving average filter. However, the flexibility of thedesign allows to use moving average filter. The median filter:

-   pads the current signal with i × B zeros;-   pads the previous signal with B zeros (this ensures equal length of    the signal whose median is calculated);-   takes the median of the resulting two signals; and-   stores the median filter as the previous signal to be used for    subsequent frame median calculation.

A resulting Median Filtered MAP estimation for a subject is shown inFIG. 20 for a length of 10 minutes. The estimated waveform has beenobtained by ConvLSTM algorithm combined with the LSTM-MODWT algorithm.

Reconstruction

The final reconstruction of the signal, as depicted in the overalldesign FIG. 13 , is the output from the LSTM-MODWT algorithm fromabove - the pulse pressure waveform, combined with (added to) the outputof the MAP, SBP or DBP waveforms. This generates the absolute bloodpressure output:

P(t) = ΔP(t) + ϕ

It will be appreciated by those persons skilled in the art that furtheraspects of the methods 100 110 will be apparent from the abovedescription of the system 10. Further, those persons skilled in the artwill also appreciate that the methods 100 110 is embodied in software,or program code, that can be supplied to the computer 14 in a number ofways, such as on a memory.

Finally, it is to be understood that various alterations, modificationsand/or additions may be introduced into the constructions andarrangements of parts previously described without departing from thespirit or ambit of the invention.

1. A method of estimating blood pressure of a subject, the methodincluding: receiving a photoplethysmogram (PPG) signal from alight-based Pulse-Plethysmography sensor applied to the skin of thesubject; processing the PPG signal using a Wavelet transform algorithmto derive PPG wavelet coefficients in a plurality of frequency bands;estimating blood pressure coefficients by processing the PPG waveletcoefficients using a machine learning algorithm that has been trained ontraining data including PPG wavelet coefficients derived from PPGs fromlight-based Pulse-Plethysmography sensors applied to the skin of testsubjects correlated with Invasive Arterial Blood Pressure (ABP) waveletcoefficients derived from simultaneously obtained ABP measurements ofsystolic and diastolic blood pressure of the test subjects; andreconstructing a pulse blood pressure waveform between systolic anddiastolic blood pressure of the subject from the blood pressurecoefficients.
 2. A method according to claim 1, further includingprocessing selected ones of the PPG wavelet coefficients for processingusing the machine learning algorithm based on an energy level present ateach of the frequency bands exceeding a threshold level.
 3. A methodaccording to claim 2, further including providing the selected ones ofthe PPG wavelet coefficients in a 2-Dimensional matrix of features tothe machine learning algorithm.
 4. A method according to claim 3,wherein each column of the 2-Dimensional matrix is autoregressiveincorporating near past and near future samples of the selected ones ofthe PPG wavelet coefficients.
 5. A method according to claim 1, whereinthe machine learning algorithm is a Recurrent Neural Network (RNN)algorithm.
 6. A method according to claim 5, wherein the RNN algorithmis a Long short-term memory (LSTM) model.
 7. A method according to claim6, further including processing the PPG wavelet coefficients usingmultiple LSTM models.
 8. A method according to claim 7, furtherincluding reconstructing output of the multiple LSTM models to form thepulse blood pressure waveform using an inverse Wavelet transformalgorithm.
 9. A method according to claim 8, wherein the Wavelettransform algorithm is a Maximally Overlapped Discrete Wavelet Transform(MODWT) algorithm, and the method further includes reconstructing thepulse blood pressure waveform using an inverse MODWT algorithm.
 10. Amethod according to claim 1, further including pre-processing the PPGsignal before processing using the machine learning algorithm bylow-pass filtering the PPG signal with a cut-off frequency for thefrequency bands.
 11. A method according to claim 1, further including:processing the pulse blood pressure waveform using a further Wavelettransform algorithm to derive pulse blood pressure wavelet coefficientsin a plurality of frequency bands; extracting features from the pulseblood pressure wavelet coefficients in the plurality of frequency bands;estimating mean arterial blood pressure (MAP), systolic blood pressure(SBP) or diastolic blood pressure (DBP) coefficients by processing thefeatures using a further machine learning algorithm that has beentrained on the training data; and reconstructing the MAP, SBP and or DBPwaveforms from the MAP, SBP and or DBP coefficients, respectively.
 12. Amethod according to claim 11, further including combining the pulseblood pressure waveform and one or more of the MAP, SBP or DBP waveformsto reconstruct an absolute blood pressure waveform of the subject.
 13. Amethod according to claim 11, wherein the further machine learningalgorithm is a Convolutional Long short-term memory (ConvLSTM) networkand the features are a 2-Dimensional feature matrix including temporaland spatial features of the pulse blood pressure wavelet coefficients.14. A method according to claim 13, further including estimatingintermediate MAP, SBP or DBP coefficients using the ConvLSTM network.15. A method according to claim 14, further including: selectingfeatures of the PPG signal and the pulse blood pressure waveform; andestimating further intermediate MAP, SBP or DBP coefficients byprocessing the features of the PPG signal and the pulse blood pressurewaveform using deep neural networks, respectively, that have beentrained on the training data.
 16. A method according to claim 15,further including: concatenating the intermediate and furtherintermediate MAP, SBP or DBP coefficients to form a concatenated outputvector; and passing the output vector through multiple layers of hiddenneurons in an output network to generate the estimated MAP, SBP or DBPcoefficients.
 17. A method according to claim 1, wherein the light-basedPulse-Plethysmography sensor includes a light source having an emissionwavelength and a photodiode having a detection wavelength, wherein theemission wavelength and the detection wavelength are around anisosbestic wavelength where oxygenated and deoxygenated blood absorbssame amount of light.
 18. A method according to claim 17, wherein thelight source is an Infrared (IR) light source.
 19. (canceled) 20.(canceled)
 21. A method of monitoring blood pressure of a subject, themethod including: applying a light-based Plethysmography sensor to theskin of the subject for a period of time; and estimating the bloodpressure of the subject according to the method of claim 1 at designatedintervals over the period time.
 22. A system for estimating bloodpressure of a subject, the system including: a light-basedPulse-Plethysmography sensor configured to be applied to the skin of thesubject to generate a photoplethysmogram (PPG) signal; a processor indata communication with the sensor; a memory; and software resident inthe memory and accessible to the processor, the software including aseries of instructions executable by the processor to configure theprocessor to: process the PPG signal using a Wavelet transform algorithmto derive PPG wavelet coefficients in a plurality of frequency bands;estimate blood pressure coefficients by processing the PPG waveletcoefficients using a machine learning algorithm that has been trained ontraining data including PPG wavelet coefficients derived from PPGs fromlight-based Pulse-Plethysmography sensors applied to the skin of testsubjects correlated with Invasive Arterial Blood Pressure (ABP) waveletcoefficients derived from simultaneously obtained ABP measurements ofsystolic and diastolic blood pressure of the test subjects; reconstructa pulse blood pressure waveform between systolic and diastolic bloodpressure of the subject from the blood pressure coefficients; and outputthe pulse blood pressure waveform to a display.
 23. (canceled) 24.(canceled)